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ABSTRACT 


Steady solutions about a slender sharp-edged delta wing in a supersonic 
freestream for moderate and high angles of attack are obtained numerically by 
time integration of the unsteady compressible three-dimensional laminar 
Navier-Stokes equations. The main features of the flow, including primary and 
secondary separation, and vortex position and strength, are adequately simu- 
lated in the numerical solutions. Improved resolution of the computational 
grid in the leading-edge region from a previous solution had considerable 
effect on the accuracy of the solutions. Good agreement between numerical 
solutions and experimental data was obtained for two cases. A local 
timestepping procedure is used to speed convergence by approximately a factor 
of two. 
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internal specific energy 
total specific energy 
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model length, mm 
freestream Mach number 
pressure 

Prandtl number, 0.73 for air 
components of heat flux vector 
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Cartesian coordinates in axial, normal, and spanwise directions, 
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numerical damping coefficient 
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finite difference step size 
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]( ' , molecular viscosity coefficient 

e , C transformed body fitted coordinates 

P density 

T xx’ T xy , T xz , 

t ,t , x components of viscous stress tensor 
yy yz zz 

Subscripts 

adiabatic wall value 
evaluated at wing leading edge 
maximum value 
minimum value 
local value 
pitot conditions 
freestream value 


INTRODUCTION 


years, with the advent of the supercomputer, computational 
aerodynamicists have devoted considerable effort to the solution of increas- 
ingly complex three-dimensional separated flows with a view toward the solution 
of airframe design problems. One of the most interesting and most practical 
of these is the separated flow field associated with a slender sharp-edged 
delta wing at angle of attack. Characterized by a pair of leeward spiraling 
vortices emanating from the separated leading-edge flow this inviscid-viscous 
flow problem has been computed using several techniques. Higher order panel 
methods (1-3) have enjoyed some measure of success by modeling the vortex 
roll-up as a sheet of quadratic doublets. However, this method requires the 
presence of a vortex and some awareness of its size and position. These 
restrictions are removed when a 3-D Euler approach is used (4,5). The physical 
meaning of such solutions is questionable, however, given the apparent 
dependency on grid resolution (6) and on the use of numerical viscosity. Using 
a conical flow approximation, Vigneron (7) et al. obtained a Navier-Stokes 
solution for a delta wing having a subsonic leading edge. With this approach, 
the location and basic structure of the primary vortex may be accurately 
computed. However, in cases where the adverse pressure gradient in the 
outboard region of the wing, induced by the leading- edge expansion, is of 
sufficient strength to produce a secondary separation, the conical approach 
becomes inadequate as evidenced by the results of Reference 7. Previous 
results presented by Rizzetta and Shang (8) illustrated that by using the full 
3-D compressible laminar Navier-Stokes equations secondary effects may be 
accurately resolved numerically. The current investigation adopts this 
technique for an additional high angle-of-attack case and attempts to enhance 
the results of the previous study by improving the grid resolution and 
leading-edge definition. The present study investigates numerically the 75° 
sweep delta wing studied experimentally by Monnerie and Werle (9) depicted in 
Figure 1. 
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GOVERNING EQUATIONS 


The governing equations for this computation are taken to be the unsteady 
compressible three-dimensional Navier-Stokes equations in mass-averaged 

variables which may be expressed notationally in the following chain-rule 
conservative form: 


au 

3t 




,u n iiv 
*3x* 3y* 3z' 


JP 




H 




3F n 

35 




8n 




3; 

3G 

35 

+f llL 

+( 3x* 

3n 

3y* 

3tu 

3z' 

3G 

3n 

+ (— 
'3x' 

35 

3y* 

11) 

3z ; 

3G 

35 

3H 



3H 




3H 

35 




3U 




35 


0 ( 1 ) 


Mere x,y,z *re Cartesian coordinates in the axial, normal and spanwise 
directions respectively and e,ri,5 are the corresponding transformed coordi- 
nates in a body-oriented system and 
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The Cartesian -velocity components are given by u, v and w, p is the density , 
E the total energy per unit mass, and e the specific internal energy. 
Completion of this system is provided by the perfect gas laws 

( 13 ) 

P = pRT 

and the Sutherland law for the molecular viscosity coefficient y. Freestream 
values are specified as boundary conditions for all dependent variables at the 
upstream, outboard, upper and lower computational boundaries , i.e. ,,U « U^. At 
the wing midspan a plane of symmetry is imposed which results in; 
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On the upper and lower wing surface 
U “ V “ v - o 
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Because only computations corresponding to supersonic freestream condi- 
tions were considered, no formal mathematical downstream boundary conditions 
were needed. However, due to the requirements of the numerical algorithm 
used, values of the dependent variables at the downstream computational 
boundary are obtained by second-order extrapolation from the interior domain. 
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NUMERICAL PROCEDURE 


The laminar computation performed in Reference 8 for M=1.95 and a=10° 
was repeated in the present study. Since it was of interest to compute a high 

angle of attacks a 20° case was chosen for the same Mach number. This 

necessitated increasing the size of the computational domain due to the 
expansion of the physical region of interest. The grid size was increased 

from 20 x 64 x 40 to 21 x 70 x 50 in the directions and was algebraically 

generated in a fashion similar (Figure 2) with two notable exceptions. In the 
work of Reference 8 the plane of symmetry was straddled to obtain a 
second-order centered difference implementation of the symmetry condition. 
While this technique was successful, the grid in this region was rather 
coarse. For the present study points were exponentially clustered near the 
plane of symmetry (Figure 2) to improve resolution in that region. This 

becomes particularly critical for the 20° case because the gradients near the 
plane of symmetry are greater due to the growth of the primary vortex and its 
inboard migration. The second modification to the grid structure was made at 
the wing tip, where resolution was increased and hence it was expected that 
the expansion about the leading edge might be more accurately resolved. This 
improvement, however, was gained at the expense of decreasing the time step 
(Table 1) due to the decreased step size in the y direction. In order to 
alleviate some of this degradation of efficiency, the lower surface n lines 
were bent down and away from the surface near the tip, thereby increasing the 
size of Ay^^ (Figure 3) . 

For the purpose of improved vortex resolution the current grid system 
employs an additional 6 q-grid lines (yielding a total of 46) on and above the 
wing surface and an additional 5 £-grid lines (total - 29) inboard of the wing 
leading edge. Further, 5 £-grid lines were added to the outboard field 
region to resolve the expansion about the leading edge. An additional 
n-C plane was added at X = 0.05 to improve the resolution near the apex. 

L 


This grid, although sized for the a=20° case, was used for both the a=10° 
and 20° cases to determine if for the present a=10° case a higher degree of 
accuracy might be obtained than its predecessor and for the purpose of 
uniformity between the current solutions. 


Steady-state solutions to Equation 1 were generated using the time- 
dependent explicit unsplit two-step predictor-corrector finite-difference 
algorithm due to MacCormack (10) which has evolved as a reliable technique for 
numerical solution of a wide variety of fluid dynamic problems. As part of 
this algorithm, a commonly employed fourth-order pressure damping term (11) 
was used to suppress numerical oscillations in regions with large gradients in 
the dependent variables. The form of this term is given as 
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which is then added to the new value of U at each time -step. In the present 
study, as in the previous effort, a damping coefficient of 3 - 2.0 was used. 
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Freestream values were used as the initial conditions save for the 
interior boundary where surface conditions were applied. In order to remove 
the initial transients, the numerical flow fields were allowed to evolve for 
100 time steps in the time accurate mode with At chosen such that the maximum 
Courant Fredicks Lewy (CFL) number was 0,5. Subsequently, the CFL number was 
increased to 0.8 and the computation was allowed to progress to steady state 
in a local timestepping fashion, where each grid point advances at its own 
rate (At-^ corresponding to a CFL number of 0.8) for a given iteration. This 
procedure has the advantage of advancing the solution to a steady-state value 
much more expeditiously since the more coarse regions of the physical domain 
are able to advance more quickly than in a time accurate mode. This modifi- 
cation resulted in the appearance of the primary and secondary vortex struc- 
tures in the numerical solution much earlier than in the work of reference 8. 

The computations were made on a Cray 1-S computer using a vectorized 
solver (12) written specifically for the Cray 1 computer. The solver advances 
the solution from one time step to the next in planes by marching in the 
([-direction, thus minimizing the data flow to and from central memory. Using 
the previously defined computational grid, approximately 1.124 x 10^ decimal 
words of storage were required. A data processing rate of 6.9 x 10~ i CPU 
seconds per time step per grid point was achieved where vectorization occurred 
in the ^-direction, with a vector length of 50. 

RESULTS 


Numerical solutions were generated for the flow about the delta wing 
illustrated in Figure 1 at a Mach number of 1.95 (0.5 normal to the leading 
edge) for 10° and 20° angles of attack. The freestream conditions correspond 
to those of the experimental work of Monnerie and Werle9, the specifics of 
which are given in the table on Figure 1. A comparison of pitot pressure 
contours from reference 8 to the current a = 10° solution is given in Figure 
4. Both solutions are a marked improvement over the results of reference 7 
which employed a conical approximation and did not capture the secondary 
separation. However, while the shape and location of the primary and 
secondary vortices are roughly the same, the pitot pressure levels on the wing 
upper surface are not. The solution of reference 8 shows less evidence of 
pitot pressure loss in the primary vortex region as shown by the large dis- 
parity in the Prp/P-p = 0.7 pressure levels. This may be attributed to the 

improved resolution 00 in the leading-edge tip region which leads to a more 
accurate computation of the leading-edge expansion, hence the improved corre- 
lation with experiment as shown in Figure 5. The expanded grid size in the n 
and ^-directions may also account for some of this improvement as the outer 
boundaries are much farther from the wing. In addition, the increased 
resolution of the apex region by the addition of the X = 0.05 station improves 

L 

the modeling of the flow being entrained into the inviscid core. 
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The improved resolution in the p or y direction accounts for the in- 
creased clarity of the upper surface "oil flow" patterns shown in Figure 6. 
While the position of the primary and secondary lines of separation (SI and 
respectively) and reattachment (A^, A ? ) are nearly the same, their sharpness, 
particularly near the trailing edge, is improved. This improvement is attrib- 
uted to a second-order treatment of the downstream boundary condition. However, 
in spite of these differences, the aerodynamic loading perceived on the wing 
is similar as shown in the Cp distributions given in Figure 7. Here the suction 
peaks associated with the primary and secondary vortices are clearly in evidence 
and are similarly located. 

The development of the cross-plane velocity with X is illustrated in 
Figures 8 and 9 for 10° and 20° respectively. Here the growth of the primary 
and secondary vortices and their inboard migration with angle of attack are 
evident and the non-conical nature of the flow field near the apex is re- 
vealed. Interesting to note is the fact that the secondary vortex originates 
at approximately X = 0.3 for both cases. These same features are also evident 
L 

in the upper surface "oil flow" pattern given in Figure 10. While the flow is 
non-conical in the apex region, analysis of Figures 7-10 indicates that for 
the majority of the planform the flow is fairly conical. 

The effect of angle of attack on loading is emphasized in the comparison 
of pressure coefficient depicted in Figure 11. The expected increase in 
loading and the inboard shift of vortex-associated surface pressure loss with 
a is evident. From analysis of the pitot pressure contours in Figure 12 one 
can observe the growth of the secondary vortex and its displacement effect on 
the primary. The increased circumferential velocity in the primary vortex 
with angle of attack induces the observed pressure loss in the vortex core. 
Comparison of the 20° solution to the experiment of Monnerie and Werle in 
Figure 13 illustrates the remarkable accuracy of the laminar Navier-Stokes 
computation in spite of the rather high Reynolds number of the experiment. 


CONCLUSIONS AND DISCUSSIONS 

Steady laminar solutions of the full 3-dimensional Navier-Stokes 
equations have been generated for the supersonic (M = 1.95) flow over a 75° 
sweep sharp-edged delta wing at 10° and 20° angles of attack. Comparison with 
experiment and with the Navier-Stokes solution of reference 8 shows the 
following: 

1. For both cases, the flow field solution correlates well with experi- 
ment. The position and strength of the primary and secondary vortices are ac- 
curately predicted. This leads one to believe that a laminar Navier-Stokes 
solution to the leading-edge vortex problem is adequate for obtaining pitot 
pressure, particularly for configurations with a Reynolds number less than 1 x 
10^ such as in the present study, at least until an adequate turbulence model 
is developed for vortical flows. 
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2. An accurate geometrical modeling of the leading edge and considerable 
resolution of the tip region appear critical to the accurate computation of 
the leading-edge expansion. The agreement obtained in the present study is 
attributed in part to these two modifications to the grid of reference 8. A 
more accurate treatment of the downstream boundary condition also accounts 
for the improved solution in the trailing-edge region. In addition, the 
expansion of the computational domain in the q and £ directions probably re- 
sulted in some improvement in the 10° solution. 

3. Because secondary separation develops axially as well as radially, 
methods employing a conical assumption will not adequately resolve a large 
portion of the flow field. In addition, resolution of the apex region is 
critical due to its highly non-conical structure. 

4. The use of a local timestepping procedure increased the convergence 
rate by roughly a factor of two. It is estimated that each solution of the 
laminar 3-dimensional compressible Navier-Stokes equations presented herein 
required approximately 2 hours of CPU time on a CRAY 1-S computer to reach 
steady state based on a four order of magnitude reduction in L2N0RM for all 
dependent variables. 
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Figure 6. Upper surface oil flow pattern, a = 10 . 
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Figure 12. Effect of a on Pitot pressure, X/L = 0.8. 
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